Correcting errors in a quantum gate with pushed ions via optimal control 
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We analyze in detail the so-called "pushing gate" for trapped ions, introducing a time dependent 
harmonic approximation for the external motion. We show how to extract the average fidelity for the 
gate from the resulting semi-classical simulations. We characterize and quantify precisely all types 
of errors coming from the quantum dynamics and reveal for the first time that slight nonlinearities 
in the ion-pushing force can have a dramatic effect on the adiabaticity of gate operation. By means 
of quantum optimal control techniques, we show how to suppress each of the resulting gate errors 
in order to reach a high fidelity compatible with scalable fault-tolerant quantum computing. 
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I. INTRODUCTION 



Trapped ultracold ions have represented a major can- 
didate for the implementation of scalable quantum in- 
formation processing since the beginning of this research 
field. The first proposal of an ion-based quantum com- 
puter by Cirac and Zoller in 1995 |l| has been followed 
by a great variety of other schemes based on ions [3L on 
other quantum optical systems like neutral atoms [1, H| 
and on solid-state systems [j| as well. With the progress 
of experimental techniques and the demonstration of en- 
tangling quantum gates based on several different candi- 
date physical systems, the focus has progressively shifted 
toward the fulfilment of scalability desiderata [f|, that is, 
the realization of quantum gates with very high fidelities, 
in the range 0.999 - 0.9999. 

Gate errors in a real implementation of a given quan- 
tum gate scheme can be reduced by different means. 
Some errors arise from (or are increased by) experimental 
imprecisions of a technical nature and can be controlled 
by careful alignment, stabilization etc. of the experimen- 
tal apparatus. Other errors stem from unaviodable in- 
teractions with the environment and can be reduced by 
simply completing the gate in as short a time as possi- 
ble. Typically, a gate scheme can be made faster by sim- 
ple scaling to higher intensities, shorter distances etc. If 
such simple optimizations of the gate prove insufficient, 
one needs to consider changes to the scheme itself and 
trade simplicity for improved performance. This is ex- 
actly the goal of quantum optimal control techniques Q 
which allow for a precise tailoring of the system's evo- 
lution by time-dependent tuning of some external pa- 
rameters. With sufficient control over these parameters, 
a given target state can often be reached with minimal 
errors even over short gate operation times. The appli- 
cation of these methods to quantum information systems 
requires in turn a very accurate simulation of the dy- 



namics and a careful understanding of the targeted error 
sources. This is precisely the aim of this paper, in the 
specific case of the two-qubit ion gate first proposed in 
and subsequently analyzed in H, H|. In this "push- 
ing gate" the qubits are encoded in the internal states of 
two ions. Each ion is held in a separate microtrap and 
state-selective "push" potentials are applied in order to 
modify the distance and thus the Coulomb interaction 
between ions, see Sec. Qll] below. We shall first point 
out a series of issues that arise when the assumption of 
spatial homogeneity of the ion-driving force is dropped, 
and subsequently develop a way to correct each of these 
issues, exploting a range of ideas, including in a crucial 
way optimal control methods. 

It should be noted from the outset that in the present 
paper, we analyze the pushing gate without what is called 
the 7r-pulse method in Rcfs. y§, Q, where it was shown 
to dramatically reduce some types of errors. The 7r-pulse 
method is a spin-echo technique and requires the gate to 
be repeated with the internal state of the ions flipped. 
Typically, single particle operations like flipping the in- 
ternal states can be done with high fidelity, and it is 
reasonable to expect that eventually the 7r-pulse method 
will be used. However, internal state control is at least in 
principle a separate issue from the pushing gate opera- 
tion itself. Keeping the design process modular in spirit, 
it is relevant to optimize the gate without this additional 
trick and thus pave an alternative road to high fidelities. 
As will become clear below, we extract quite general noise 
reduction methods from the automated numerical opti- 
mization results and it is an interesting topic for future 
research to combine these with the 7r-pulse method. 

The remainder of the paper is organized as follows. In 
Sec. |TT] we introduce the general setting of conditional 
dynamics gate schemes, a useful approximation for sim- 
ulating such a gate, and a measure of the gate errors. 
In Sec. Mil we specialize to the pushing gate. The un- 
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optimized performance of the gate is reported in Sec. IIVI 
and in Sec.[V]we show how this performance can be signif- 
icantly improved by a combination of "manual" changes 
and numerical optimal control methods. Finally we con- 
clude in Sec. I VII The Appendix contains a number of 
more technical results and derivations. 



II. CONDITIONAL DYNAMICS 

The basic idea of the so-called pushing gate is one of 
conditional dynamics, i.e. we apply potentials depend- 
ing on the internal state of the two ions. The inter- 
nal states themselves are not changed during the gate, 
or, in the case of the 7r-pulse method, are changed on 
a much shorter time-scale than the external dynamics. 
This means that the analysis of the problem splits into 
four separate evolutions for the external state, one for 
each of the logical (internal) states, 00, 01, 10, and 11. 
In the following, these four evolutions will be denoted 
"branches" and will be indexed by /3 € {00,01,10,11}. 
The complete Hamiltonian can be written in the form of 
a sum of internal^external factorized terms 



H tot (t) 



E 



Hf{t) 



(1) 



and it results in an evolution operator of a similar form 
U*°\t) = J2\PM®W(t). (2) 



Ideally, when t = T at the end of the gate, the Up x 
should differ from each other by at most a phase factor 
multiplying a common unitary operator U° 
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Uj?(T) =e ie ?U[ 
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so that U tot itself can be factorized: 




(3) 



(4) 



The internal evolution is then that of a phase gate, while 
the external evolution can in principle be undone using 
internal-state independent potentials. 

The requirement ([3]) is very hard to achieve and would 
make the gate completely independent on the initial ex- 
ternal state. However, we can typically assume to have 
some degree of control over the initial external state, e.g. 
by cooling the particles before the gate. This means 
that need only hold when restricted to a subset of the 
complete Hilbert space, typically the states of relatively 
low energy. In Appendix [B] we show how to evaluate the 
performance of the gate in general. For now, we note 
that since only the low energy part of £7J X (T) will be im- 
portant, we can focus on getting a good approximation 
for this part when trying to simulate the gate dynamics. 



A low initial energy means particles localized near the 
potential minimum and this suggests using a harmonic 
approximation to the real potential [22}. The simplest 
choice is to Taylor expand around a fixed point, which is 
not changed during the gate operation. The next level of 
refinement is to expand around the instantaneous poten- 
tial minimum. This works very well if the gate operation 
is nearly adiabatic so the particles stay near the (moving) 
minimum at all times. However, it may be desirable to 
make fast and substantial changes to the potential during 
the gate and that may induce pronounced non-adiabatic 
dynamics. In that case, the harmonic approximation can 
still be a good one provided it is done around the classical 
trajectories of the particles. Typically these trajectories 
cannot be computed analytically, but for any moderate 
number of particles it is a numerically simple task to find 
them. Below we will use this method and show how it 
leads to a relatively simple characterization of E/J x . 



A. Harmonic approximation 

In this section we focus on a single branch of the evolu- 
tion and thus suppress the (3 index. Let us denote by x(t) 
the classical trajectory, which is found by solving classi- 
cal equations of motion. The time-dependent, second 
order Taylor expansion of the potential V(t, x) around 
x(i) reads simply 

7 so (*,x) - F (0) (<)+Ax T V«(<) + iAx T y( 2 )(i)Ax, (5 ) 



with Ax = [x — x(i)] and 



VM(t) = V(t,x), 



(0) 



d 2 V 

dxidxj 



(t,x). 



Note that we will still use x as our coordinate, i.e. we 
are not changing to a coordinate system moving with 
x(f), we simply use a potential that approximates the 
real potential close to x(t). Collecting terms of equal 
order in x leads to the alternative form 

V so (t, x) - E(t) - x T F(t) + ix T X(t)x, (7) 

with 

E(t) = V<°)(t) -x T (i)V«(t) + ix T (i)W 2 )(i)x(t), 

F(i) = -vW(i) + y( 2 '(«)x(t), 

K{t) = V {2) {t). 

(8) 
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B. Gaussian evolution 



C. Fidelity 



The big advantage of choosing a second order approx- 
imation to the real potential is that this restricts the 
corresponding approximate /7 ex (i) to be Gaussian for all 
t. Let us introduce the compact notation q = (x, p) and 
define the matrix J by 



J = 



In 
-In 



(9) 



where n is the number of degrees of freedom. Then the 
usual canonical commutation relations can be written as 



(10) 



With the potential of Eq. ([7]) and a matrix of parti- 
cle masses M — diag(mi, . . . ,m n ) the time-dependent 
Hamiltonian becomes 



1 



H so = -p L M- L p + V so (t,it). 



(11) 



In Appendix [X] we show that such a Hamiltonian leads 
to an evolution operator of the form 



U(t) 



(12) 



where D c is a displacement operator and Wj, a squeezing 
operator 



-ic Jq 



-i^q 6q 



(13) 



The scalar (f>, the vector c, and the matrix S — exp(J6) 
should satisfy the following equations of motion: 







1 



d_ 

dt 



c = Jhc 



(14) 



dt 



JhS, 



where the 2n x 2n matrix h is defined by 
hit) \ K ® ° " 



(15) 



The form of solution ([T2"j) ~ (fl"4")) holds for any second 
order Hamiltonian. In the particular case where V so is a 
Taylor expansion of a real potential around the classical 
trajectory x(i), the equation of motion for c reduces to 
the exact equation of motion for q = (x, p) where p is the 
classical momentum. In the following, we will therefore 
write q instead of c. It is then important to remember 
that the right-hand sides of Eqs. (fl4l are in general non- 
linear functions of x. 



We can quantify the performance of the gate by calcu- 
lating the average fidelity -F avg (see Appendix[B]) between 
the obtained output state and the ideal one when the in- 
put state is varied. One can then separate out three kinds 
of contributions to the deviation of f avg from 1 (the per- 
fect gate) 



1 



Favg — Eg 



Eq + Es (16) 

The three types of errors each have their physical inter- 
pretation. The most straightforward one pertains to the 
sloshing errors Eq which correspond to a residual mo- 
tion of the ions after the gate has been completed and 
the micro-traps are again at rest. The phase errors Eg 
are errors in the gate phase. Finally, the breathing er- 
rors Eg are induced by differences in the harmonic ap- 
proximation parameters around the classical trajectory 
for different internal states. For example, in the case we 
consider, when the particles are pushed closer together, 
the second order term in the Coulomb repulsion becomes 
larger, cf. Eq. (|23|) below. 

In our model we assume that systematic, local phase 
errors can be undone. Then an explicit calculation in 
Appendix |Bl shows that the phase errors are given by 

Eg = ^(9 m -6 al -9 w + 9 11 -Trf (17) 

with 9p = —(f>p+ Tr[b/3j], where 7 is the covariance matrix 
of the external state, see Section IB 31 of the Appendix. 
Note the inclusion of the Tr[6 J g7] terms in the definition of 
9p: These terms correspond to phase contributions from 
average excitation energy in the traps and are therefore 
temperature dependent through 7. For our parameters 
they are small. 

The sloshing errors are given by 



E q = 20 51 ^ ~ q /3) T ^( c l" _ ^p) 



(18) 



Q</3 



and the breathing errors are 



~ £ Tr [(b a - bf,)7(b a - bp)j] 



1 

160 



Lj2Tr[(b a -b p )J(b a -b p )J\. 



(19) 



a</3 



Here 0/3, q», and bp are defined in Eqs. (|12I13|) . and are 
the variables describing the Gaussian approximation in 
the j3 branch of the evolution. Again 7 introduce tem- 
perature dependence. 



III. THE PUSHING GATE 

Let us now focus on the particular case of the pushing 
gate. Here we have two ions, each in a separate micro- 
trap [23l | . To further simplify the discussion, we concen- 
trate on just one spatial dimension, i.e. there are two 
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degrees of freedom, n = 2. The ions are assumed to be of 
identical mass, mi = m 2 = m, and thus M = mS.2- The 
potential energy consists of a micro-trap for each ion, 
time- and internal state-dependent pushing potentials, 
and the Coulomb interaction. The pushing potentials 
can be realised as optical dipole potentials generated by 
focused laser beams. The time-dependence of these po- 
tentials is most easily achieved by controlling the inten- 
sity of the laser and the state-selectivity by polarization 
selection rules [9]. We assume that the form for the in- 
ternal state labeled by /3 is 



mw 2 ij + 



i=l,2 



Attcq \d + x 2 — xi\ 



E fru) 



i=l,2 



a 



(20) 



The trapping potentials are assumed to be perfectly har- 
monic. The state-dependent pushing amplitudes are 
such that the ions are only pushed if they are in the in- 
ternal state 1, f^\t) = 5p iA f(t). Note that the 
factor means that ion 1 is pushed to the left and ion 2 
to the right. Non- linear contributions to the pushing po- 
tentials are included via the constant G. The harmonic 
oscillator ground state size is do = y/h/mcu. The two 
coordinates x\ and xi are taken to have origin in the 
respective trap centers, a distance d apart. Experimen- 
tally, theparameters in (|20l) can be varied quite a lot 
(see e.g. [9]). Trap distances d from lOO^m all the way 
down to 1/xm are within technological reach. The trap 
frequency u> can be chosen in the range 2ir x 10 4 — 10 7 Hz 
which for e.g. Ca + ions will mean an oscillator length a 
from 150nm down to 5nm. 

It should be noted that the use of optical dipole po- 
tentials to generate the state selective pushing forces will 
in general introduce large single qubit phases due to ac- 
Stark shifts. This is not a problem as such, but it means 
that even small fluctuations in laser intensities will lead 
to loss of gate fidelity. For the particular case of the 
pushing gate, the ac-Stark shifts can be balanced against 
the Coulomb energy as discussed in Ref. [9J. Obviously 
ac-Stark shifts are common to many gate proposals that 
uses optical potentials. An experimentally demanding, 
but quite general solution is to compensate the shifts 
along the lines of Refs. [13, [HI- In the present work, we 
focus on errors that are more directly related to the mo- 
tion of the ions and assume that the push potentials are 
effectively non- fluctuating. 



A. Dimensionless Hamiltonian 

The relative strength of the Coulomb interaction to the 
trapping potentials turns out to be conveniently quanti- 



fied by 



a 7r<L d 

j 2 d 2 d? fuv 



(21) 



which is the ratio of the energy scale of the Coulomb and 
trap potential energies at the equilibrium positions of 
the ions. In Ref. it was found that e -C 1 is the most 
promising regime. In oscillator units, the Hamiltonian 
for the branch labeled by f3 reads 



1 



■£ 5 

i=l,2 



4|l + ^(i 2 -x 1 )| 



(22) 



B. Harmonic approximation 

When Eqs. © are specialized to the pushing gate, we 
get the following: 




(23) 



if>(t) = -(-l)Vf>( i)T i^ 

4 ao 



1 + 



4«(<) = i+2/r(*)G 



1 

+ 2 £ 



1 



_o ( 

d \ 



1 



1 



d \ 



x m - x m 



(24) 
(25) 

(26) 



When these expressions are inserted into Eqs. ([T4]) and 
(p~5]) we are ready to simulate the gate. 



IV. RESULTS OF SIMULATION 

A. Choice of parameters 

Even with the simplifications we have introduced, 
there are still a lot of parameters in the problem. The 
optimal "working point" will always be dependent on ex- 
perimental considerations beyond the simplified model 
treated here. For a discussion of parameters and design 
decisions, see Ref. @. For concreteness we have cho- 
sen to focus on a limited set of parameters. We first of 
all assume the individual ion traps to be very well sep- 
arated and let a/d = 0.001 in all calculations. Likewise, 
we assume a reasonably low value for e of 0.04. Such 
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parameters would result from e.g. Ca ions placed in 
micro-traps with trapping frequencies of uj ~ 27r x 5MHz 
and separated by a distance of ~ 7^m. For a the trav- 
eling wave configuration with beam waist w considered 
in Ref. [9], G = 4:(a/w)(w/2xo — 2xq/w) where xo is the 
initial position of the ion relative to the beam center. 
For realistic focusing of the push beam, w ~ l/im this 
suggests to vary the non-linearity coefficient G between 
and 3 x 1 CP 2 . As the initial temporal shape of the push 
pulse we choose a Gaussian f(t) = £exp(— t 2 /r 2 ), where 
the amplitude £ should be chosen to give a gate phase 
of 7r. A simple estimate (for G = 0) suggests that we 
choose [1] 

= ; 7 = ■ (27) 

v /V8e v /l + e/2wT 

The temporal width of the pulse, r, should be within an 
order of magnitude from the trap period if we want a fast 
gate. We will mainly look at r in the range 1 to 10 trap 
periods, which for the parameters quoted above results 
a maximum excursion due to the push in the range from 
12ao down to 3ao- 

B. Phase errors 




FIG. 1: Phase errors with non-uniform pushing forces. The 
physical parameters are: e = 0.04, a/d = 0.001, and uit = 3.5 
(black), 5.5 (red) and 7.5 (blue). Results for temperatures of 
T — 0.125, 1 and 8 x hu/ks, are plotted for each value of r 
and are indistinguishable from each other. 



The choice of push amplitude expressed by Eq. (|27|) is 
not optimal. This can be seen in Fig. [T] where we plot 
Eg as a function of G. Even for G = the gate phase 
is not exactly it. For lot — 3.5 we see that a nonzero 
G can improve the gate phase. This is not surprizing, 
but also not very useful as we shall see below that Eg 
is in general easy to reduce. In Fig. [1] results for three 
different temperatures are plotted, but the dependence 
on temperature is completely negligible. 

C. Sloshing errors 

Let us now turn to the errors described by the Eq term 
in Eq. (|C8|) . the sloshing errors. Fig. [5] shows how these 
errors are strongly dependent on G, the strength of the 
non-linearity of the pushing potential. A series of minima 
of Eq- as a function G can be seen. The optimal values of 
G depends on the chosen duration of the pulse, r. Each 
minima is associated with the ions performing an integer 
number of "non-adiabatic oscillations" during the push 
pulse. This is illustrated in Fig. [3] where the trajectory 
of ion 1 with respect to its trap minimum is plotted for 
values of G that are below, at and above the one that 
leads to the lowest Eq. 

In contrast to Eg above, Eq depends noticeably on 
whether T is 0.125, 1 or 8 x hu/k-Q- Higher temperatures 
always increase the sloshing errors and for k-^T '/hu 1 
we find that Eq scales as T since 7 does [see Eqs. (|B21[) 
and (lB22|) ]. 

Curves for three different values of r are plotted in 
Fig. [5] and it is immediately clear that one can dramat- 




0.005 0.01 0.015 0.02 0.025 0.03 

G 



FIG. 2: Sloshing errors with non-uniform pushing forces. The 
physical parameters are: e = 0.04, a/d = 0.001, and uit = 3.5 
(black), 5.5 (red), and 7.5 (blue). We plot Eq as a function 
of G. As explained in the text, there exist nonzero values 
of G where the sloshing is strongly suppressed. Results for 
temperatures of T = 0.125, 1 and 8 x foj/fe are plotted as 
dashed, full and dotted lines, respectively. 
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0.3 




-0.4 -0.2 0.2 



^cl,l *^min,l 

FIG. 3: Phase-space trajectories for ion 1 for a push duration 
of ojt = 5.5. The coordinates are relative to the potential min- 
imum in which ion 1 is trapped. Perfect adiabatic evolution 
would correspond to the ion simply following this minimum 
and thus the trajectory would be the single point (0, 0) in this 
plot. Since the push is not infinitely slow, the ion will first lag 
behind the moving minimum and then oscillate in the moving 
potential. When the potential minimum again approaches its 
original position, the ion may happen to have just the right 
position and speed in order to end up at rest. Whether this 
is the case depends (for fixed push pulse) on G: the higher 
G, the more the confinement is increased during the push. In 
the figure, the thick line corresponds to G — 0.002, which is 
nearly optimal w.r.t. returning the ion to rest. The thin lines 
correspond to G = 0,0.001,0.003,0.004. The lowest G gives 
the outermost curve over the main part of the "loop" in the 
figure. 

ically decrease sloshing errors by make the gate slower 
and thus more adiabatic. The suppression of Eq- is expo- 
nential and this is thus in general an efficient strategy. 



D. Breathing errors 

We now turn to the Eg term of Eq. (fTB)) . These 
"breathing" errors come from the different changes in the 
effective quadratic Hamiltonian for the different branches 
of the evolution. In Fig. |4] we plot Es for different values 
of G and r and for different temperatures of the exter- 
nal motion. Results for temperatures of T =0.125, 1 and 
8xfiw/fcB are shown and it is first of all clear that Eg 
depends strongly on T. We also see that Eg is nearly 
proportional to G 2 t. That larger G leads to larger errors 
is not surprising, but that larger r does is rather counter- 
intuitive: larger r means a more smooth and thus more 
adiabatic push. It also means a smaller amplitude for the 
push since the ions will have more time to pick up the 
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0.002 0.004 0.006 0.008 0.01 

G 

FIG. 4: Breathing errors with non-uniform pushing forces. 
The physical parameters are: e — 0.04, a/d — 0.001, and cur — 
3.5 (black), 5.5 (red), and 7.5 (blue). Results for temperatures 
of T = 0.125, 1 and 8 x hw/k^ are plotted as dashed, full and 
dotted lines, respectively. For all three temperatures, larger 
r leads to larger Es. Es is increasing approximately linearly 
with G 2 t, i.e. for fixed G, Es will be proportional to r, the 
pulse duration! This is very different from the behavior of Eq 
above, which is rapidly decreased by increasing r and thereby 
making the push more adiabatic. 



gate phase, cf. Eq. (f2~T|) . Let us discuss the explanation 
for this behaviour in more detail. 

For exponential suppression of errors to be valid, the 
evolution should be well into the adiabatic regime. At 
first sight, the relevant timescale is the oscilla- 

tion period of the micro traps. Since the traps are 
assumed to be far apart (a <C d) the parameter e is 
small and the normal modes of the system have peri- 
ods shifted little from this value: the CM mode is in 
fact unaffected by the Coulomb interaction and has fre- 
quency oj while the relative motion mode oscillates at 
VI + e oj. With lot 3> 1 one should therefore not be able 
to put excitations into either of these modes. However, 
it is perfectly possible to transfer excitations between 
the modes, as the adiabatic timescale for this process is 
(VI + e w — w ) 1 ~ e~ 1 aj~ 1 /2. Such a transfer will be 
induced by mixing of the CM and relative motion during 
the gate operation. A linear push potential will not mix 
the two, but a non-linear one will. 

Another effect to remember is that when the instanta- 
neous oscillation frequencies change during the push, the 
external motion will pick up different phases depending 
on the number of excitation quanta. Like the transfer of 
excitations, this effect of course disappears if the system 
is cooled to the ground state. A perturbative calculation 
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to lowest order in a/d, G and e gives the result 



{ 1 



sinri (fuv/2k B T) 
1 

sinh 2 (fiw/2fc B T) 



1 + exp 



exp(-2arVH. (28) 



The prefactor gives the scaling behaviour both in the 
"naive" nonadiabatic limit lot < 1 and in the more rele- 
vant intermediate region 1 < ujt < e _1 : 



E s cx £ z G 2 uj 2 t 2 cx 



G 2 ujt 



(29) 



In fact, even in the adiabatic limit ewr 3> 1 this scal- 
ing holds true since one term in Eq. (|28p does not con- 
tain an exponential damping factor with r 2 . This un- 
suppressed term is stemming from the above mentioned 
effect of time-varying instantaneous mode frequencies. 

From Eq. (|28p we can also understand the strong 
temperature dependence of the breathing errors. For 
k-QT/hw 3> 1, the breathing errors will scale approxi- 
mately like T 2 . However, as seen in Figure 01 high tem- 
peratures require very low values for G. For low tem- 
peratures, note that one term in Eq. (|28|) is not sup- 
pressed even at T — 0. This term stems from changes in 
the ground state widths of the two instantaneous normal 
modes and is adiabatically suppressed when ujt 1. To 
find the dominant term for very low temperatures and 
short pulses one should do a higher order perturbativc 
calculation. 



V. OPTIMIZING GATE PERFORMANCE 



From the simulations in Sec. UVI we learn that without 
improvement high fidelities require either very low values 
of G or cooling the external motion almost to the ground 
state. In this section we shall see how a better perfor- 
mance can be achieved by modifying the temporal shape 
of the push-pulses. 




0.008 0.01 



FIG. 5: Improving Ee by iteratively adjusting the push- 
amplitude f. The pulse-duration is r = 7.5 and the other 
parameters are as in the previous figures. Four sets of curves 
are shown, each with results for T = 0.125 (dashed), 1 (full), 
and 8xhu)/kB (dotted). The uppermost set (looks like a single 
curve) is for the un-optimized pulse-amplitude and identical 
to the r = 7.5 curves of Fig. [T] The progressively lower sets 
are for one and two iterations of the amplitude scaling de- 
scribed in the text. Note that as the error gets smaller, tem- 
perature begins to have an effect. This simple adjustment 
is capable of reducing Eg below 10 -6 for all the considered 
G < 0.01. 



plying this algorithm to the t — 7.5 curves of Fig. [T] As 
can be seen, Eg is rapidly reduced and can be brought be- 
low e.g. 10 -6 in a very modest number of iterations. For 
simplicity we ignore the temperature dependent Tr[&7] 
contribution to Eg when rescaling the pulse. This is the 
reason for the k^T = 8Huj (dotted curves) departing from 
the k^T — 0.l25huj and k^T — IHuj curves especially at 
low G. 



B. Fast gate: eliminating sloshing in x 



A. Correcting the phase 

Our first step will be to correct the gate phase by a 
simple scaling of the push pulse-shape. From Fig. [1] we 
know that typical errors can be well above the percent 
level. Our strategy is based on the observation that the 
simplest estimate of the gate-phase suggest that it scales 
as the square of the push-amplitude, £. [A general gate- 
phase replaces the tt in the numerator of Eq. (l27|) .] We 
therefore divide £ by the square-root of the ratio of the 
observed gate-phase and the ideal gate-phase (tt) and re- 
peat the propagation. In Fig. [S] we show the results of ap- 



A big advantage of the simple harmonic approximation 
is that it becomes feasible to solve the equations of mo- 
tion many times with different temporal shapes of f(t) 
in order to optimize the performance of the gate. Rather 
than simple trial-and-error we will apply the global con- 
trol algorithm of Krotov, which is guar anteed to improve 
the performance at each iteration |12h14| . The relevant 
equations for our case are given in Appendix [C] 

In general, it is desirable to complete the gate in as 
short a time as possible. This will limit many undesired 
effects and it will ultimately enable faster quantum com- 
putations. A fast gate, however, means that the pushing 
force will deliver a rather abrupt impulse. This can lead 
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iterations 



FIG. 6: Optimization of pulse-shape to eliminate sloshing. 
We plot Eq, the contribution of sloshing to the total infidelity, 
as a function of the number of Krotov iterations performed. 
The parameters in this example are: e = 0.04, a/d = 0.001, 
G = 2 x 10" 5 . Each curve corresponds to a separate value 
of the pulse duration, lot =3.5, 5.5, 7.5, with larger r always 
giving a lower value of Eq. 



to excitations of the external motion being left after the 
completion of the gate, limiting the fidelity. In this sec- 
tion we show how such "sloshing" effects can be avoided 
by using optimal control. 

We start from an initial Gaussian temporal shape of 
the push. The overall amplitude is first optimized itera- 
tively to get the desired gate phase as described above. 
We then run the Krotov algorithm to get a better shape 
of the pulse. We assume a non-uniform pushing force, 
G = 2 x 10~ 3 . The result is plotted in Fig. H As can be 
seen, the influence of sloshing motion can be decreased 
by a couple of orders of magnitude in a modest number 
of iterations. 

To investigate the physical mechanism behind the re- 
duction of the sloshing error, we plot in the lower panel of 
Fig. [7] the difference between the optimized pulse f op t(t) 
and the original Gaussian pulse /o(i) = £ exp(— t 2 /r 2 ). 
This difference looks a lot like a simple cosine-wave 
with a period close to 2ttlo~ 1 multiplied by a Gaus- 
sian of the same width as fo. Thus the optimized 
pulse is approximately of the form f op t(t) ~ fo(t) + 
Acos(uit) exp(— t 2 /r 2 ). An approximative calculation of 
the sloshing excitation for the simplest G = case re- 
veals that for such a pulse, the non-resonant contribu- 
tion of the bare Gaussian pulse is cancelled by a resonant 
contribution from the cosine-modulated pulse. Since the 
resonant response is much stronger, only a small, nega- 
tive A is needed for this cancellation. More precisely, the 
optimal A from first order pertubation theory is given 
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FIG. 7: Upper panel: pulse shape both before (full line) and 
after (dashed line) Krotov optimization. Curves for three dif- 
ferent pulse duration are shown, lot — 3.5, 5.5 and 7.5, other 
parameters are as in Fig. [5] Only for the shortest pulse is the 
optimized curve distinguishable from the original Gaussian. 
Lower panel: Difference between optimized pulse and initial 
Gaussian pulse (after adjusting overall amplitude to reduce 
phase errors). The left plot is for G = 0, the right plot is 
for G = 2 x f0 -3 . Each curve has been normalized to the 
prediction of a G = perturbative calculation, which is seen 
to describe well the G = case as all curves have a maximal 
excursion of approximately -f, while the G / case is only 
qualitatively similar. 

by — i£ exp(— lo 2 t 2 /A) and Fig. [7] shows that this is also 
what the Krotov algorithm converges to for G = 0. The 
"strategy" of the Krotov algorithm in this case seems 
therefore to be well understood. For G ^ 0, it is more 
difficult to predict the value of A, but nonetheless the 
Krotov algorithm seems to be highly efficient. 

C. Minimizing breathing errors 

We now know that phase errors and sloshing errors 
can be controlled and we turn to the breathing errors 
of Fig. SJ Without optimization, these errors put rather 
stringent limits on the parameters. In order to keep E$ at 
an acceptable level, either very low temperature or very 
small G is required. For very low temperature k^T — 
0A25Hlo, we need just G < 10~ 2 to get E s below 10~ 4 , 
but if we assume a more modest cooling to k^T = ltko 
the same error-level require G < 2 x 10~ 4 . Note that even 
at G = 0, breathing errors persist and that for k^T — 
8hio they never get below the 10 -5 level. These errors 
stem from the high order terms in the Coulomb potential 
which have been ignored in Eq. (l28l) . 

As described in Sec. IIVD[ breathing errors cannot be 
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eliminated by simply increasing the push duration r: 
First of all the adiabatic time-scale is ~ e _1 a; /2 which 
will mean a slow gate and secondly even in that limit er- 
rors from the change in normal modes frequencies remain 
and even increase, cf. discussion below Eq. (l29|) . It turns 
out that a simple application of the Krotov algorithm is 
also not very efficient in reducing the breathing errors. A 
partial explanation for this can be found from the pertur- 
bative calculation leading to Eq. (f28|) and the analysis of 
sloshing error-reduction above: since the adiabatic time- 
scale for the breathing errors is long, the Gaussian pulses 
we consider are not adiabatic w.r.t. breathing errors and 
thus the admixture of a small resonant component in the 
push-pulse will not be enough to get the cancellation we 
found in the case of sloshing errors. In fact, the amplitude 
of the cosine modulation should be comparable to the to- 
tal amplitude for the relatively short pulses considered. 
There is nothing to be gained from a small amplitude 
modulation and thus the "linear" version of the Krotov 
algorithm we apply (see Appendix [C]) will not work. 

In fact, in order to cancel out the contribution to the 
breathing errors due to mode frequency changes, sign 
changes in the push amplitudes are required during the 
pulse. This is beyond the simplest physical implementa- 
tions where the push amplitude is proportional to some 
laser intensity. In principle it is possible to play with de- 
tunings to implement the sign changes, and this will in 
fact give many of the advantages of the 7r-pulse method, 
see Ref. @. Allowing negative push amplitudes and 
putting "by hand" an optimized cos(ewt/2) modulated 
contribution, we have been able to e.g. reduce Eg below 
1(T 6 for ujt = 7.5, G = 2 x 1(T 3 , and T = 1 x huj/k B - 
Compared to the results reported in Fig. SI this is a re- 
duction by more than 3 orders of magnitude. Unfortu- 
nately, the strongly modified pulse now gives rise to large 
sloshing errors. To obtain an overall satisfactory fidelity 
we use a combined strategy: we first put the breathing 
error reduction by hand, then iteratively reduce phase 
errors and finally use the Krotov algorithm to reduce the 
sloshing errors. In Fig. [5] we show results of this strategy 
starting from lot — 7.5. 



VI. CONCLUSION 

In this paper we have shown how a time-dependent, 
quadratic approximation to the Hamiltonian can be a 
useful tool when analyzing quantum gates based on con- 
ditional external dynamics. The resulting equations 
are much more manageable than the original two-body 
Schrodinger equations. This is especially true if one were 
to include more spatial dimensions than the one consid- 
ered here: A full time-dependent, three dimensional, two- 
body wavefunction calculation is an extremely demand- 
ing numerical task, whereas the corresponding quadratic 
approximation will be much more manageable. 

We used the developed method to show how to im- 
prove on a naive design of the pushing gate. This was 



io- 2 


X 

V 


QgK8® 




io" 4 






X 


io" 6 




V 






o 


□□□□ 


□□□□□□□□□□□□□□□□ 


io" 8 




v v 


V 







6 


Krotov 



FIG. 8: Combined strategy for reducing infidelity for e = 0.04 
and G = 2 x 10~ 3 . The three contributions Eg, E q , and Es 
are labeled by V, o, and □, respectively. Their sum is la- 
beled by x . The leftmost column are the results for a simple 
Gaussian with ujt — 7.5. The breathing errors dominate. In 
the next column, breathing errors have been reduced (but 
sloshing increased) by using a cosine-modulated pulse based 
on Eq. (|28p . The overall amplitude is then iteratively opti- 
mized to reduce phase errors as described in Sec. IV Al As 
can be seen, 3 iterations are more than sufficient to render Eg 
completely insignificant. The dominating error type is now 
sloshing and in the third column, the Krotov algorithm is 
applied to finally reduced the total infidelity below 10 -4 . 



done including a non-uniform contribution to the pushing 
force. An important lesson of our analysis and simula- 
tions, is to pay attention to changes in the symmetry 
of the Hamiltonian during the gate operation. In the 
present case, a nonlinear push potential invalidates the 
separation of the dynamics into CM and relative motion. 
This opens up another type of non-adiabaticity, namely 
transfer of excitations between the two normal modes. 
The adiabaticity parameter for this type of error is clot 
and for small e, adiabaticity will require gate times much 
larger than the charateristic time of the micro-traps. An 
efficient counter-measure is to decrease temperature so 
that there are in fact no excitations to transfer between 
modes. Failing that, one should increase e as much as 
possible and, somewhat counter-intuitively, do the gate 
as fast possible. Sloshing errors puts a lower limit on 
the gate-time, but as we show an optimized choice of the 
temporal shape of the push can dramatically reduce this 
problem. 

By analyzing the way that the Krotov optimized pulse 
reduces sloshing errors, we identified the basic mechanism 
as a destructive interference between the non- resonant, 
non-adiabatic contribution from the finite push pulse du- 
ration and a resonant contribution from a small ampli- 
tude superposed oscillation of the push force. Generaliz- 
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ing this idea to deal with breathing errors, we were able 
to reduce them by several orders of magnitude. However, 
since breathing errors are not significantly adiabatically 
suppressed for the considered pulse durations, the de- 
structive interference required sign-changes in the push 
amplitude, which introduce experimental complications. 
The suppression of breathing errors also came at the price 
of increased sloshing errors, but we showed that the Kro- 
tov algorithm once again was able to improve the pulse 
shape. 

One may ask to what extent is the final pulse shape 
optimal for the given overall gate-time? This is an in- 
teresting question in general and in this study we saw 
examples both where the Krotov algorithm seemed to 
exhaust the potential in its "strategy" (the G = case 
in Fig. [7]) and where it was not able to find an optimiza- 
tion. In the latter case we could improve the pulse by 
hand (eliminate the breathing errors by destructive in- 
terference). The problem of optimality is related to the 
question of a quantum speed limit (QSLj [l5| and this 
connection has been studied in Ref. Note, how- 

ever, that in our case the hamiltonian is time-dependent 
and what we want is in fact to leave the external motion 
unaffected after the pulse. It would be interesting to in- 
vestigate such a general adiabaticity problem along the 
same lines as the work on the QSL. 
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Appendix A: Evolution under quadratic Hamiltonian 

In this section we show that a second order Hamilto- 
nian leads to an evolution operator that can be written 
like U(t) Eq. (fT2")l . For alternative parameterizations, see 
Refs.[13, [HI- We will do a direct calculation showing 
that U(t) fulfills the Schrodinger equation 



for the displacement operator and 



d - 

i—U = H so U. 
at 



(Al) 



The demanding part of the the calculation involves dif- 
ferentiating expontials of time dependent operators. A 
useful formular can be found in e.g. [l9| and involves in- 
tegration over an auxiliary variable r\. It results in 

i^- t £>c=J o A,c c T Jq £»t c d V b c 



c T J (q — 77c) drj D c 



1 



(A2) 



=c 1 J(q--c]D c 
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ifWb = J W nb -q T 6qH^ dr) W b 

=~ J (e-i Jb qfb (e- r ' Jb q) dr, W b 

=h T J T f e^Jbe^Jb dv qw b 
2 Jo 

-54 TjT (!«")•- 4* 



(A3) 



for the squeezing operator. It is now easy to show that 
the equations of motion (fT4"]l for <fr, c and S — exp(Jb) 
leads to U = exp(-i^)D c W b fulfilling Eq. (JXTJ [HJ. 



Appendix B: Fidelity for Gaussian evolutions 

In this section we derive expressions for the fidelity 
as a function of the variables used to characterize the 
evolution in the harmonic approximation, tf>p, q~g, and 
S , £ = 00,01,10,11. 

If we assume that the initial state of the system is a 
product of an internal state density matrix and an exter- 
nal state density matrix, p®<r, we get the following state 
for the internal degrees of freedom after the application 
of U tot of Eq. ©: 



pi =Tr c 



U tot p ( 



a/3 



(Bl) 



Here "o" denotes the element-wise matrix product (the 
Hadamard product) and R is the matrix given by: 



[R] af s = Tr 



U, 



> (c>r) 1 



(B2) 



It is easy to see that R er Hermitian and that all its 
diagonal elements are 1. In particular, Tri? = 4. Slightly 
less obvious is it that R is positive semi-definite: Let 
c e C 4 . Then: 



c' Rc = Tr 
= Tr 



[{£^}H£c;i>°> 



(B3) 



where we have used the cyclic property of the trace and 
the fact that the trace of a product of positive semi- 
definite operators is non-negative. 

The element-wise product form in Eq. IB 1 1) is per- 
haps not the most illuminating. If we diagonalize R = 
J2h w hwj l , we get instead a Krauss operator sum form: 



h=l 



(B4) 
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with K h = diag(wh). 

There are different ways to define the fidelity of the 
gate. The one used in Refs. [|| and [|| is as the minimum 
fidelity of the obtained final state w.r.t. the wanted final 
state when the input state is varied. This means that 



^ mi „ = min(^|C/ t / 9;C/o|V') 



(B5) 



with = IV'XV'I an d Uq the gate operator we aim for. 
In the case of a phase gate, Uo is diagonal in the logical 
state basis: 



U 



E 





and we get the simpler minimization problem: 

F min = miiip T Rp 
{pi} 



(B6) 



(B7) 



where thep € and Y^Pi = The matrix R = UoRUq 
is Hermitian, but since p is confined to be real, only its 
real symmetric part contributes. The minimization in 
Eq. (|B7p is a so-called quadratic programming problem 
and very efficient numerical methods for its solution ex- 
ist. Given R it is therefore simple to calculate F m - m on 
the computer. However, a more direct evaluation is pos- 
sible if fidelity is instead defined as an average over input 
states as 



iv 



WOj^oW dV. 



(B8) 



Here S denotes the normalized states (unit sphere) in 
C™ and the volume element dV is such that J g2 „_i dV — 
1. For -Fkvg a compact formula exist [20] and using it in 
the present case leads to 



Favg ~ 4(4 



(B9) 



and AD = D-(D) 



where (D) = Tr [DU™ m a cx (UZj 
Calculating R to second order in the AD's, we get: 



R a p =exp {iA9 a ) 



x 1 



x exp (— iA6p) 



ADj - 2AD a AD 



(B12) 



This form is useful, as it separates the infidelity into sys- 
tematic phase-errors (the exp(±iA0 a ) factors) and "de- 
coherence" (factor in curly brackets). The phase errors 
A9p = (Dp) — Qp can be made small by tuning the av- 
erage of laser powers etc. and this can usually be done 
very well. The challenge will therefore most often be to 
suppress the fluctuations, i.e., the terms in curly brackets 
in Eq. (|B12|) . 

Assuming also the A#'s to be small, Eq. (|B10j) be- 
comes: 



AD 



l-^avg = ^ £ (W - A6 P ] 2 + /(AD, 

(B13) 

At first sight, this form might seem dubious since only 
differences in the A#'s and AI?'s enter. However, one 
should remember that any common evolution on the four 
branches can be absorbed into U^ m in Eq. (|BTT|) . This 
emphasizes that for the implementation of a single gate 
on the logical state, the external motion must not neces- 
sarily be returned to its initial state as long as the final 
state is common to all logical input states. Typically, 
the further requirement that energy is not pumped into 
the external degrees of freedom by repeated application 
of the gate must be made. In the particular case of the 
pushing gate, this requirement is in fact already hidden in 
Eq. (IB13I) since the (3 = 00 branch contains no pushing. 
In other cases, one could apply cooling to the external 
state between gate operations. 



or when using the properties of R: 



2. The non-local part of the phase 



1 



l-F avg = -^(l-Re R a0 



a</3 



(BIO) 



1. General small errors 

Typically, we will be mostly interested in situations 
where the four logical states leads to almost identical 
evolutions for the external states. It is then useful to 
write 



Uf=e^{iD p )U{ 



cx 
com 



: exp 



(i{bp)) exp(iADp)U° 



(Bll) 



We are seeking to implement the phase- gate (|B6|) . In 
many cases, the Op 's are not so important individually 
since single-particle operations are easy to perform and 
only the truely non-local phase is interesting. Assuming 
that perfect single-particle phase changes can be imple- 
mented on average, it is straightforward to show that one 
should replace ^2 a< a[A9 a — AOp] 2 by 



[A0oo - a# 01 - A6 10 + Ae n y 



(B14) 



in Eq. (IBT3)) . 

This simplified view of single-particle phase-changes 
should of course be revisited in a more complete anal- 
ysis of any given proposal for quantum-computing. In 
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the present work we use the replacement (|B14[) through- 
out, but let us emphasize that fluctuations in the single- 
particle phaserotations are more naturally incorporated 
in the AD terms of Eq. (|B13|) than in the AO terms: One 
simply model the fluctations as a consequence of some 
fluctuating parameter which can be included in <7 ex - 



Appendix C: The Krotov Algorithm 

Optimizing the temporal shape of the push pulse is 
done using the Krotov algorithm 1121 . For an introduc- 
tion to the method, see e.g. Ref. |13j . 



The Gaussian case 



1. Auxiliary variables 



For Gaussian evolutions like (fT2j) and a Gaussian (e.g. 
thermal) external state a with covariance matrix 

Hi =2 + mi) - (li) (&) / B15 ) 
=ReTr fe^cr] - Tr [q t a] Tr [qja] 

and vanishing means 

(g<) = Tr [ qi a] = (B16) 

one gets phase contributions 

(Dp) = -<f>p - Tr [bp-y] (B17) 

and decoherence terms 

((A£> Q - ADpf) = (c a - c ) T J 7 J T (c Q - cp) 

+ ^Tr [(b a - bp)~/(b a - 6,9)7] 

+ ^Tr [(b a - bp)J(b a - bp) J] (B18) 
o 

In general, for a harmonic oscillator in thermal equi- 
librium at temperature T, the covariance matrix is given 
by 



1 



1 



Vthermal 



2 tanh 



2k B T 



_A_ 

fvmuj 



(B19) 



where ks is the Boltzmann constant. In the present case, 
the CM and the relative motion are separately in ther- 
mal equilibrium and for the corresponding dimensionless 
position and momentum operators [xqm = (%i + ^2 )/2 
etc.], we get: 



The key ingredient in this approach is a function 
$(t, {4>p,qp, Sp} p=oa,oi,iQ,ii) which allow us to translate 
the global goal of improving the final C7 tot (T) to a local 
problem of choosing a better f(t) for each t. Construct- 
ing $ is in general very difficult, but it is relatively simple 
to get a linear approximation to it. The coefficients in 
this approximation will constitute a set of auxiliary vari- 
ables. For each branch, the equations of motion for the 
auxiliary variables <f>p, qp and Sp are determined by the 
requirement that they are conjugate to the physical vari- 
ables 4>p, q^, and Sp, respectively. 

Let us focus on a single branch and suppres the (3 
index like in Sec. Ill Al We then need to construct 
T-L(t] <p, q, S; 4>, q, S) such that Eqs. (fT4|) can be written 



d 

d_ 



d 



Is- 
dt 13 



dq~i 

d 



n. 



(CI) 
(C2) 
(C3) 



This leads simply to 



ft(t;0,q, S; & q,S) 



1 



£(t,x)--F^,x) 



-T 

q 



J/i(t,x)q- JF(t,x) 
Tr S T Jh(t,x.)S 



(C4) 



Then the equations of motion for the auxiliary variables 
become 



with 



7rel 



and 



^ tanh - ^77^ Tri 

2k B l 



1 

7CM = - 



7 = 7CM © 7rel 







1 



2 tanh 



2k B T 



\ 

2 



(B20) 



(B21) 



(B22) 



In the limit e <C 1, a we have approximately 7 oc I4 if we 
use the set of individual-ion operators (xx, X2,P\,Pi)- 



d ~ d 
d 

^S = VsU = h(t,ii)JS. 



(C5) 
(C6) 
(C7) 



The equation of motion for q is rather involved since H 
depends on x in a complicated manner through h, F and 
E. It can be rewritten as two coupled time-dependent, 
forced harmonic oscillators. Note on the other hand that 
(j) is time-independent and that JS solves the same equa- 
tion as S. 
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2. Objective function 



terminal conditions, i.e. boundary conditions at t = T: 



Our ultimate goal is to improve the fidelity of the gate. 
However, it is somewhat impractical to apply this as the 
objective in the Krotov algorithm: Calculating the fi- 
delity is only simple for small errors and in general it 
depends on e.g. the temperature of the external motion. 
Instead we shall work with a simpler function of the vari- 
ables 0, q and S for the four branches. The reduction of 
this objective function should tend to increase the fidelity 
of the gate. Based on the fact that in the pushing gate, 
the branch (3 = 00 is not subject to any time-dependent 
forces, we choose the following: 

J ({<A8) q,3, 5*/3}^=00, 01,10,ll) = J<j> + Jq + Js 

1 2 



'oo — Wl 



§5>( 



?10 + 011 - 7T 
2 



X/3 - X 00 



Pp - Poo 



(C8) 



~>I3 — <J00 



)> 



(10 



The term with 0's aim to ensure the correct phase in the 
phasegate, while the other terms aim at identical evolu- 
tion for the external motin in the four branches. In the 
limit e -C 1 Eq. (IB13|) formally justifies the use of our 
chosen objective function, given the extra proviso that 
we are only interested in the non-local part of the phase. 



= -(-l) 73 (0oo - 0oi - 0io + 0ii - t) | T 

MT) = -Vx,j| T 

[3x o - xqi - xio - In ) , ,5 = 00 

J T 



MT) = - 



Sp{T) = - 



3Poo - P01 ~ P10 ~ P11) 
Pa ~ Poo 



^s? J\ T 

3Sqq — Sqi — Sxo — on 

Sf) — Sqq 



(C9) 



(CIO) 



p = 00 
p ^ 00 



(Cll) 



p = 00 

p^oo 



(C12) 



3. Terminal conditions 

The objective function supplements the auxiliary- 
variable equations of motion (|C5HC7p with the following 



where (— l)* 3 is +1 for P = 00 and 11 and —1 for p = 
01 and 10. These equations express the values of the 
auxiliary varibles at time t — T in terms of the physical 
varibles also at t = T and give the input to the backwards 
propagation of the auxiliary variables, cf. Ref. [l2| . 
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